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Abstract 


This  paper  describes  an  efficient  Monte  Carlo  sampling  plan 
for  estimating  the  distribution  of  maximum  flow  in  a  directed 
network  whose  arcs  have  random  capacities.  Such  a  network  can  be 
used  to  represent  a  multi-state  system  whose  multi-state 
components  are  subject  to  deterioration  in  capacity  by  random 
amounts  at  random  points  in  time.  The  proposed  sampling  plan 
uses  an  easily  computed  a  priori  upper  bound  on  the  complementary 
distribution  function  to  obtain  an  unbiased  point  estimator  with 
smaller  variance  than  crude  Monte  Carlo  sampling  allows. 
The  paper  also  describes  procedures  for  interval  estimation  and 
for  assessing  when  the  sampling  experiment  has  achieved  a 
specified  accuracy.  To  facilitate  sampling,  the  paper  presents 
a  characterization  of  deterioration  based  on  cumulative 
processes,  leading  to  the  treatment  of  arc  capacities  as  being 
m u 1 1  i  n  o  r  m  a  1  1 y  distributed.  A  technique  is  described  for  checking 
the  appropriateness  of  this  model  with  regard  to  lower  and  upper 
bounds  on  capacity.  A  procedure  is  also  described  for  deriving  a 
confidence  interval  on  the  measure  used  to  assess  variance 
reduction.  An  example  illustrates  the  sampling  plan  and  a  concise 
summary  gives  all  steps  needed  to  implement  the  plan. 


Key  words:  Maximum  flow,  Monte  Carlo  sampling,  multi-state 
reliability 
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Introduction 


This  paper  describes  a  method  for  estimating  the 
complementary  distribution  function  of  maximum  flow  in  a 
multi-state  reliability  system  whose  components  have  capacities 
that  decrease  randomly  through  time  and  for  which  maximum  flow 
characterizes  system  performance.  River  systems  with 
accumulating  silt  deposits  and  distributed  manufacturing  systems 
whose  individual  processing  units  deterioriate  through  time 
provide  two  examples  of  such  multi-state  systems.  Although  a 
considerable  literature  has  developed  on  the  reliability  of 
multi-state  systems  with  multi-state  components  that  exhibit 
random  d e t e r  i  o  r  a  t  i  o n  in  performance,  this  literature  mostly 
concentrates  on  characterizing  properties  of  performance  measures 
induced  by  the  inherent  features  of  broad  classes  of  multi-state 
systems.  For  example,  see  Barlow  and  Wu  (  1  978  ),  Baxter  and  Kim 
(1984),  Block  and  Savits  (1984)  and  Mak  (1985).  A  considerably 
more  modest  literature  exists  for  the  equally  important  topic  of 
computing  multi-state  performance  measures.  For  example,  see 
Butler  (1982),  Hudson  and  Kaper  (1985)  and  Kulkarni  and  Adlakha 
(1985). 

The  present  paper  is  a  contribution  to  this  computational 
literature.  Specifically,  it  characterizes  system  performance  by 
the  concept  of  maximum  flow  and  describes  a  Monte  Carlo  sampling 
plan  for  estimating  the  complementary  distribution  of  maximum 
flow  with  considerably  better  statistical  accuracy  than  crude 
Monte  Carlo  sampling  would  allow  for  the  same  cost.  For 


convenience  of  exposition,  we  use  the  nomenclature  of  network 
analysis  to  characterize  a  system.  Although  the  emphasis 
throughout  this  paper  is  on  system  reliability,  the  methods 
described  here  clearly  apply  more  generally  to  the  considerably 
broader  class  of  maximum  flow  problems  with  capacities  subject  to 
random  variation. 

Let  G  =  (  V ,  E , s , t  )  denote  a  directed  network  with  node  set 

V,  arc  set  E  =  |l,...,n),  source  node  s  and  sink  node  t 

s.teV.  At  time  t£0,  arcs  which  represent  the  components  of  the 

system  have  capacities  B-|(t) . Bn(t).  The  capacity  of  each  arc 

is  subject  to  deterioration  so  that  Bj(t)  2  Bj(t+A)  A£0  ISiSn. 
The  times  at  which  the  deteriorations  occur  are  random,  as  are 
the  magnitudes  of  the  incremental  deteriorations  themselves.  Let 
T  denote  the  set  of  all  minimal  s-t  cutsets  of  G  and  define 
the  capacity  of  cutset  C  at  time  t  as 

Z(t  ,C)  =  l  B. (t)  CeT. 

ieC  1 

If  capacities  were  d e t e r m i n  i  s t  i  c  rather  than  random,  then  a 
measure  of  system  status  for  each  time  t  would  be 

J1(t)  *  min  Z(t,C).  (1) 

CeT 

As  is  well  known  from  the  max-flow  min-cut  theorem  of  Ford  and 
Fulkerson  (1956),  A ( t )  is  the  maximum  possible  flow  from  source 
node  s  to  sink  node  t  at  time  t.  Then  the  function  { A ( t  )  t£0} 


shows  the  pattern  of  deterioration  in  performance  in  G  as  time 
evolves. 

When  capacities  are  random,  (1)  no  longer  provides  an 
adequate  measure  of  system  performance  since  { A ( t  )  t£0}  is  now 

a  stochastic  process.  Then  for  fixed  t£0  interest  focuses  on 

L(x,t)  =  pr[A(i)>x]  x£0.  (2) 

Also,  for  fixed  x20  interest  may  focus  on  the  first  passage  time 
tx  =  min[z:  A  (  z  )  S  x  ] 

for  which 

pr  (  tx>t  )  =  L(  x  ,  -t )  .  (3) 

Then  for  each  i 10  ( L ( x  ,  t  )  ,  x20j  gives  the  complementary 

distribution  function  of  maximum  flow  and  for  each  xSO  |l(x,t), 
x20|  gives  the  first  passage  time  complementary  distribution 
function. 

Other  measures  of  system  performance  also  are  of  interest. 
For  example, 

L(x,C,t)  =  pr[Z(t,C)  =  A ( t  )  and  A(x)>x]  (4) 

denotes  the  probability  that  C  is  the  critical  minimal  s-t 
cutset  and  that  the  maximum  flow  exceeds  x  at  time  t.  Then 
L  (  C  ,  t|  x  )  =  L(  x  ,C  ,  t  )/L(  x  ,  t  )  (5) 

=  probability  that  C  is  the  critical  minimal  s-t 
cutset,  given  that  the  maximum  flow  exceeds  x  at 


time  t  . 


Relatively,  little  work  has  appeared  on  the  numerical 
evaluation  of  the  measures  (2)  through  (5).  Somers  (1982) 
describes  a  procedure  for  computing  (2)  exactly  for  the  cases  of 
one  and  two  arcs  with  random  capacities  and  n-1  and  n-2  arcs  with 
fixed  capacities,  respectively,  and  notes  that  the  computation 
time  for  n  randomly  capacitated  arcs  grows  exponentially  with  n. 
Assuming  independent  exponentially  distributed  arc  capacities  in 
a  planar  network,  Kulkarni  and  Adlakha  (1985)  show  how  to  compute 
(2)  exactly  in  0(  |v|  •  |e|  )time,  using  the  max-flow  algorithm  of 
Itai  and  Shiloach  (1979). 

Frank  and  Frisch  (1971)  provide  a  comprehensive  discussion 
of  the  randomly  capacitated  maximum  flow  problem  with  regard  to 
the  computation  of  (2).  Their  discussion,  which  is  principally  in 
terms  of  normally  distributed  arc  capacities,  gives  lower  and 
upper  bounds.  Because  of  the  difficulty  of  computation  for  even 
moderately  sized  networks,  they  suggest  a  crude  Monte  Carlo 
sampling  plan  to  estimate  (2). 

The  present  paper  describes  a  Monte  Carlo  sampling  plan  for 
estimating  (2)  that  employs  an  easily  computed  upper  bound  to 
gain  its  computational  advantage.  In  particular,  this  bound 
enables  one  to  modify  sampling  in  a  way  that  reduces  the  sampling 
variation  inherent  in  each  trial,  thereby  reducing  the  sample 
size  required  to  achieve  a  specified  accuracy  of  estimate  when 
compared  to  crude  Monte  Carlo  sampling.  Although  the  paper 
concentrates  exclusively  on  (2),  one  should  note  that  repeated 
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application  of  the  proposed  sampling  plan  for  a  selected  set  of 
times  i  would  also  allow  estimation  of  (3). 

To  set  the  stage  for  the  proposed  method,  Section  1 
describes  a  prototypical  crude  Monte  Carlo  sampling  plan  that 
remains  the  basis  for  comparison  throughout  the  paper.  Section  2 
presents  an  improved  sampling  plan  that  exploits  a  priori 
information  on  the  upper  bound  to  improve  statistical  efficiency. 
The  plan  includes  point  and  interval  estimation.  Section  3 
describes  a  particular  model  of  component  capacity  deterioration 
in  a  network  that  enables  one  to  regard  capacities  on  individul 
arcs  as  being  mul t i normal  1 y  distributed.  This  model  greatly 
facilitates  the  use  of  the  proposed  sampling  plan.  The  section 
also  describes  how  to  assess  the  appropriateness  of  the  normal 
model  with  regard  to  lower  and  upper  bounds  on  capacity. 

Section  4  presents  a  procedure  to  guide  the  choice  of  when 
to  stop  sampling  and  Section  5  describes  inferential  methods  for 
assessing  the  extent  of  variance  reduction.  An  example  in 
Section  6  illustrates  the  proposed  procedures  in  detail  and 
Section  7  summarizes  the  essential  steps  in  Sections  2  through  b 
to  implement  the  proposed  sampling  plan. 

1 .  Crude  Monte  Carlo  Sampling 

Procedure  MFCRUDE  describes  the  steps  to  be  performed  in 
carrying  out  a  crude  Monte  Carlo  sampling  plan  designed  to 
estimate  |L(x,t)  x  c  X  }  from  data  on  K  independent  trials  or 


replications.  Here  L  ^  (  x  ,  t  )  is  an  unbiased  estimate  of 
L ( x , t  )  with 

var  Lk(x,t)  =  L(x,t)  [ 1 -L(x  ,  t )  ]/K . 

Moreover,  V(Lk(x,t))  in  step  III  is  an  unbiased  estimator  of 
var  Lr ( x , t  )  .  Note  that  in  step  II  sampling  usually  takes  0(n) 
time  per  replication.  However,  the  dominant  time  consumer  is  the 
determination  of  maximum  flow;  for  example,  taking  0(  |v|  log  |v|  ) 
time  per  replication  for  a  planar  network  (Itai  and  Shiloach 
1979). 

With  regard  to  computing  the  cell  number  j  in  step  II, 

choosing  the  flow  set  as  x^  =  a  +  8i  a  ,  $  >  0  i  «1 . r  offers 

an  advantage.  Then  for  continuous  A(t),  one  has 
j  =  [_m  i  n  {  r  ,  [  A  (  t  )  -  a  ]  +  /  6  }  J  so  that  the  time  to  determine  j  is 

constant  and  independent  of  r.  If  the  set  of  flows  X  does  not 
admit  the  representation  x^  =  a  +  Bi,  it  always  will  be  beneficial 
to  augment  the  set  so  that  it  does. 

Procedure  MFCRUDE  represents  a  baseline  sampling  experiment 
against  which  one  needs  to  compare  any  alternative  candidate 
sampling  experiment  with  regard  to  the  presumably  smaller 
resulting  variances  of  interest  and  the  usually  larger  computing 
time  required  per  replication  to  realize  these  variances. 
Section  2  describes  one  such  alternative  sampling  plan  based  on 
an  easily  derivable  upper  bound  on  (l(x,t)}. 


Procedure  MFCRUDE 

Purpose :  To  estimate  (2)  for  specified  time  t  and  selected  flow 

levels  X  =  {  x0  =  0< x i < . . . < xr } . 

Input:  Network  G  =  (  V  ,  E= (  1  ,  .  ,  .  , n )  , s  t  t  )  ,  time  t,  joint 

distribution  of  capacities,  flow  set  X  and  number  of 
independent  replications  K. 

Output:  {  L  k  (  x  ,  t  )  ,  V (  L  k  (  x  ,  t  )  )  ;  xeX}  as  estimators  of 

{  L ( x  ,  t  )  ,  v  a  r  L  k ( x  ,  t  )  ;  x  e  X  }  ,  and  {  A  S  (  x  )  xeXj  as  an 
available  input  to  continued  sampling. 

Method: 


I.  Initialization 

Set  AS(Xj[)  =  0  i  =  0,  1  ,  .  .  .  ,r. 

II.  Sampling  Experiment 

On  each  of  K  independent  replications:  sample  the 
capacities  B  i  (  x  ),...,  B  n  (  t  )  from  the  joint  distribution; 
determine  the  maximum  flow  A(x);  compute  j  =  max[i: 

Xi<A(-r)  i  =  0  ,  1 . n  ]  ;  AS(xj)  =  AS(xj)  +  1. 

III.  Computation  of  Summary  Statistics 


LK( x .  ,  t  ) 


AS  (  x  ) 


1  <i<r 


V(Lk(x.,t))  =  Lk ( x  .  ,  t ) [ 1  - Lk( x  .  ,  t )  ]/ ( K- 1  )  . 


End  of  Procedure 


Exploiting  Bounds 


2  . 

The  plan  to  be  described  here  concentrates  sampling  in  a 

region  of  the  sample  space  of  B-j(x),...,Bn(x)  where  each 

collected  sample  outcome  contains  considerably  more  information 

than  in  the  crude  Monte  Carlo  case.  Let 

H  (  x , C , t  )  =  {  Z  ( C , t ) >x  } 

so  that  (2)  has  the  form 

L ( x , t  )  =  pr [ n  H ( x , C  ,  t  )  ] 

CeT 

=  M(x,t)  N(x,t) 

where  for  =  I^Cx.-t)  £  T 

M  (  x  , -t  )  =  pr  [  H  H  (  x  ,  C  ,  x  )|  Pi  H(x,C,t)] 

Cer~r1  Cer1 

and 

N  (  x  ,  t  )  =  pr  [  O  H  (  x  ,  C  ,  x  )  ] . 


Let 

H  (x  ,  x  )  =  n  H ( x  ,  C  ,  x  ) 

C  e  r  1 

and  suppose  that  N(x,x)  can  be  computed  exactly.  Then  to 
estimate  L(x,x)  one  proceeds  as  follows: 

On  each  of  K  independent  replications  sample 
B ( x )  =  ( B1 ( x ) , . . . , Bn ( x ) )  subject  to  the  restriction  H  1  ( x  ,  x  )  ; 
determine  the  maximum  flow  A(x);  set  $  =  I  .  .  ( A  (  x  )  )  where 

\  X  j  ®  / 


I  denotes  the  indicator  function  on  the  set  A. 


T 
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Then 

-  ,  .  1  r  (k) 

M  (  x  ,  T  )  =  -  l  4) 

K  K  k  =  1 

where  the  superscript  (k)  denotes  trial  k,  is  an  unbiased 
estimator  of  M(x,t)  and 

L  (  x  ,  x  )  =  N ( x  ,  t  )  M  ( x  ,  t  )  (6) 

A  A. 

is  an  unbiased  estimator  of  L(x,t)  with 

var  L  (x,t)  =  L(x,x)  [N(x,t)  -  L(x,t)]/K.  (7) 

A 

Observe  that 

var  ( x  ,  t  )  i  N ( x  ,  t  )  / 4K 


which  can  be  a  considerably  tighter  upper  bound  than  that  of  1/JJK 
for  var  Lk(x,t).  More  importantly,  note  that 

v  a  r  L  (  x  ,  -r  )  1  -  L  (  x  ,  x  )  1 

- _£ -  =  -  >  -  ,  (8) 

var  Lk(x,t)  N(x,t)  -  L(x,x)  N(x,x) 

indicating  that  to  achieve 


var  L . 


=  var  L„  ( x  ,  x  ) 
K2 


the  ratio  K i / K 2  must  be  greater  than  or  equal  to  1/N(x,x). 


The  variance  ratio  in  (8)  is  merely  one  measure  of  benefit. 
Let  T  and  T  denote  the  mean  times  required  to  perform  one 
replication  with  crude  sampling  and  with  the  bound  respectively. 
Then  T/T  provides  a  second  measure  of  performance  and  is  usually 
less  than  unity.  As  a  single  measure  of  merit  the  quantity 

-  var  L  ( x , t  ) 

R(x,t)  =  -4~  *  - ~ -  (9) 

T  var  LK( x  ,  t  ) 

is  commonly  used  in  the  Monte  Carlo  literature.  This  quantity, 
called  the  variance  reduction  ratio,  measures  the  amount  of  time 
one  would  have  to  run  the  crude  Monte  Carlo  experiment  to  achieve 
the  same  variance  as  one  would  achieve  using  the  proposed  method 
in  one  unit  of  time.  Clearly  R(x,t)>1  indicates  that  the 
proposed  method  is  preferred.  Section  5  describes  inferential 
methods  of  estimating  (9)  from  sample  data  at  the  end  of  an 
experiment  . 

Choosing  Hj (x  ,  t  ) 

In  selecting  H  i  ( x , t  )  one  needs  to  consider  several  issues. 
First,  N ( x , t )  must  be  computable.  If  arc  capacities  are 
statistically  independent  and  the  cutsets  in  Ti(x,t)  are 
edge-disjoint,  then 


N  {  x  ,  t  ) 


n  pr [H ( x , C , t ) ] 

Ce r  ( x , t ) 


( i  o) 


whose  computation  takes  0(|  r | ( x , t )  |)  time,  once  the 
pr [H ( x , C , t ) ]  are  determined. 

To  estimate  L  (  x  , -r )  for  just  one  flow  value  x,  one  may  want 
to  choose  r i  ( x , t )  to  be  the  set  of  ed ge- d i s j oi  n t  minimal  s-t 
cutsets  that  minimizes  (10).  Unfortunately,  the  determination  of 
this  r^x.t)  is  NP-complete.  Moreover,  since  one  customarily 
estimates  the  complementary  distribution  function  for  several 
flow  values,  it  becomes  necessary  on  each  replication  to  insure 
that  H  -|  (  x  ,  x  )  holds  for  each  x  in  X.  If  r^x.T)  contains  more 
than  one  cutset,  this  can  entail  considerable  resampling  thereby 
increasing  computation  time  substantially.  To  reduce  this 
sampling  frequency,  we  describe  an  alternative  method. 

Let  r*  denote  a  set  of  ed ge- d i s j o i nt  minimal  s-t  cutsets, 
let  X  =  {Xo  =  °<xl  <  *  •  .<xr }  denote  the  set  of  flows  of  interest  and 
define 

i-  1 

P .  ( D  ,  t  )  =  pr {h( x  , D  ,  t  )  n t  n  H(x  ,C( x  ,  x  )  ,  x  )  ]  }  Der* 

i  i  j  =  i  J  J 


where 

C  (  x  t  ,  t  )  =  {  C  :  CeT*  I  P^C.t)  2  P^D.x)  DeT*}.  (11) 
Here  C(x.,x)  minimizes  P^D.x)  at  flow  x^  Now  redefine 

i 

N ( x ,  , t )  -  pr[  n  H ( x . , C ( x , , t ) , t  )  ]  (12) 

j  =  1  J  J 


which  is  again  an  upper  bound  on  L(x  ,t).  While  it  is  true  that 


(13) 


-  1  2- 

N  (  x  ,  i  )  >  pr  [  Pi  H  (  x  , C  ,  t  )  ]  , 

1  Cer*  1 

using  { C ( x  i  ,  t  )  1  S i S  r  )  ,  as  defined  in  (11),  makes  resampling 
necessary  for  at  most  one  constraint  H  (xj  ,C ( xj , i )  , t )  at  step  i, 
in  contrast  to  the  potential  sampling  for  the  |r*|  constraints 

n  H(x  ,  C , x )  that  are  active  on  each  step  if  one  chooses  to  use 
Cer* 

all  cutsets  in  r  *  .  Note  that  one  can  determine  a  set  of 

edge-disjoint  cutsets  r  *  =  |  C  . . C  j  }  in  0(n)  time  where  J  is 

the  size  of  smallest  minimal  s-t  path.  Procedure  A  shows  how  to 
determine  |c(xj,i)  1£iSr|  given  f  *  . 

Procedure  MF BOUNDS  shows  how  to  perform  K  independent 
replications  to  estimate  (2)  for  |c(xi,x),  N(x1,t)1  as  determined 
in  Procedure  A.  Note  that  at  step  i  and  i  constraints 
i 

H  H ( x  , C ( x  , t ) , t )  are  in  force,  although  only  meeting 
j  =  1  J  J 

H ( x  j  , C ( x  j  , t )  , t )  may  induce  resampling.  Since  the  probability  of 
resampling  at  step  i£2  is 

pr [ Z( C ( x . , x ) , t ) >x  ] 

1  pr  L  Z(C ( xi , t ) , t  )  >x j  (  .  j  ]  ’ 

where 

j(i)  =  max  {k:  C(x  ,x)  =  C(x.,x)}  if  C  (  x  ,  x  )  e ( C ( x  ,  i  )  l<k<i| 

1  <  k  <  i  K  1  IK 


0 


otherwise, 
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the  mean  number  of  sampled  capacities  on  each  trial  in  Procedure 
MFBOUND  is 


r 

n  +  l  C  (  x  ,  x  ) 
i  =2 


pr  [  Z  ( C  (  x  .  ,  x ) , t ) >x  .  ] 


pr[Z(C(x.,T),T)>x^ 


rp| 


and  an  upper  bound  on  the  mean  number  of  maximum  flow 
determinations  is 


r  pr[Z(C(x,T),t)>x] 
y  - i - i - 

l  :2  pr [ Z(C ( x.  , T ) , t ) >Xj  (  ^  ] 


By  contrast,  Procedure  MFCRUDE  samples  n  capacities  and 
determines  maximum  flow  just  once  per  replication.  Therefore, 
Procedure  MFBOUND  S  clearly  takes  more  time  per  replication  than 
Procedure  MFCRUDE  does. 


Procedure  A 


Purpose:  To  compute  an  upper  bound  for  {pr[A(x)>x]  xeX(. 

Input:  Edge- di s j oi nt  minimal  s-t  cutsets  r*  =  {C],...,Cj}, 

flow  set  X  =  { Xq  =  0< x i < . . . < xr  }  ,  time  x  and 
{ pr [H ( x i ,Cj ,x)  ];  ISiSr,  1 S j  < J  j  - 
Output :  {c(xifx),  NCxj.x)  1£i£r}. 

Method:  Find  the  smallest  m  that  minimizes  pr  [H  ( x -j  ,  Cm  ,  x  )  ]  for 

all  ISmSJ;  set  C(xpt)  =  Cm;  set  N(x-|,x)  = 
pr[H(x1,C(xi,x),x];  i=1. 

Until  i =r : 

A.  Find  the  smallest  m  that  minimizes 

i 

pr)H(x  ,C  ,T)n[  n  H(X  ,C(x  for  all  ISmSJ; 

1  I  Ul  j  _  i  J  J 

set  C(x.+1 ,x)  =  Cm; 

set  N(xi  +  1,x)  =  N ( x  i  ,  x  )  pr[H(x  . +1  ,C(x.  +  1  ,x  )  ,x)  ]. 

B.  If  Cm  e  (C(xj,x)  1  £  j  £  i  }  :  determine  the  largest  j 
such  that  C(xj,x)  =  Cm;  set 

N(xi+1,x)  =  N(Xi  +  i .x  )/pr[H(xj  ,C(xj  ,x)  ,x)  ]. 

C.  i  =  i  +  1. 


End  of  procedure 


Purpose:  To  estimate  (  2 )  for  specified  time  t  and  selected  _ _ J 

- ' -  , 

.  *  -  -  „  f 

flow  levels  X  =  |  x0  =  0<  x  ?  <  .  .  .  <  xr  }  . 

Input :  Network  G  =  (  V  ,  E  =  jl,...,n)1s,t),  time  t,  joint 


distribution  of  capacities,  flow  set  X, 
|C(x,t),N(x,t)  as  determined  in  Procedure  A;  xeX) 
and  number  of  independent  replications  K. 

Output :  { L K ( x , x )  ,  V ( L K ( x , x ) ) ;  x  e  X  }  as  estimators  of  (l(x,t), 

var  L  ( x , t  )  ;  xeX  }  and  ls(x)  xeX  x*x.},  as  an  available 
K  u 

input  to  continued  sampl i ng. 

Method : 

I.  Initialization 

Set  C(xqit)  =  C(x1(t);  set  S(xj)  =  0  1Si<r. 

II.  Sampling  Experiment 

For  each  of  K  independent  replications  do: 

A.  Sample  capacities  j Bj { t  )  1Si£n)  from  the  joint  distribu¬ 
tion  subject  to  the  restriction  H ( x 1  ,  C ( xq  ,  t  )  ,  t  )  . 

B.  Compute  maximum  flow  A(t). 

C.  Set  Z  =  Z(t,C(x0,t)  )  . 

D.  For  i  =  1 . r  do: 

1.6=0. 

2.  If  A ( x ) >X|  ,  6  =  1. 


3.  Otherwise: 


a.  If  C(x.,t)  *  C ( x i _ 1  , t ) ,  set  Z  =  Z ( t  ,C ( xt ,  T  )  )  . 

b.  If  Z  £  x.:  resample  (b^Ct)  jeC(xi>i)}  subject  to 

H(  x  .  ,C  (  x  .  t  )  ,  t  )  ;  set  Z  =  l  B  .  ( t ) ;  compute 

1  1  jeC(x.,T)  J 

maximum  flow  A(t);  if  A(t)>x^  6=1. 

4.  S ( x . )  =  S(x.)  ♦  6 . 

III.  Compute  summary  statistics 

L„  (  x . , x )  =  N ( x .  , t )  S(x  )/K 
&  J  J  J 


1£j<r  . 


V(L„(x,,t))  =  N2(x.,t)  [S( x  .  )/K  ]  C 1  - S ( x  .  )/K]/(K-1  )  . 
“  J  J  J  J 


End  of  Procedure 

Confidence  Intervals 

While  variances  give  an  indication  of  statistical  accuracy 
some  studies  call  for  the  construction  of  confidence  intervals 
for  estimates  of  interest  at  a  specified  level  of  significance 
a  0<a<1.  Let  L  =  L(x,t),  N  =  N(x,t)  and  S  =  S(x)  for  xeX 
where  S(»)  is  given  in  Procedure  MFBOUNDS.  Since  S(x)  has  the 
binomial  distribution  with  parameters  K  and  L(x,t)/N(x,t),  it  is 
in  principle  possible  to  derive  exact  confidence  intervals  for 
L(x,t).  However,  as  pointed  out  in  detail  in  Fishman  (1984),  the 


numerical  evaluation  of  the  interval  becomes  increasingly 


difficult  as  K  increases.  Here  we  describe  two  alternatives. 

For  binomial  samples,  the  results  in  Okamoto  (1958)  enable 
one  to  show  that  ( N  ■]  ,  N1P2)  covers  L  with  probability  exceeding  a 
where  \p  ^  <  1P2  are  the  solutions  of  the  equation 

(  K-  S )  1  n  (  1  -  4/ )  +S  in*  »  ln[(1-a)/2]  +  (K-S)ln(l-S/K)+S  ln(S/K).  (Ha) 
For  a  >  .9  these  bounds  are  considerably  tighter  than  those  that 
Chebyshev's  inequality  produces.  Also,  the  interval 


NS  +  N32/2  ±  BN  [  B2/ 4  +  S(K-S)/K]1/z 
K  +  82 


(14b) 


asymptotically  (K-*®)  covers  L  with  probability  a  and 


6  =  |  y : 


-1/ 

(2n)  2 


_x2/2 
e  dx 


-  ( 1 +o)/2  ) 


(15) 


Employing  the  tighter  normal  result  requires  care.  If 
L(x,i)/N(x,t)  is  close  to  zero  or  unity  then  convergence  to 
normality  is  slow  and  the  normal  interval  can  be  misleading.  One 
indicator  of  the  extent  of  convergence  is  the  skewness  measure 

E[L(x,t)-L(x,t)]^  N(x,t)-2L(x,t) 

o  =  - - -  =  - rr  (16) 

[var  Lk(x,t)]3/2  {  K  L  (  x  ,  t  )  [  N  (  x  ,  t  )  -  L  (  x  ,  t  )  ]  (  2 

A  skewness  close  to  zero  supports  the  contention  of  normality. 
Section  6  uses  a  sample  value  for  Q  to  make  this  assessment. 


Occasionally  it  may  be  of  interest  to  compute  simultaneous 
confidence  intervals  for,  say,  (l(x,t)  xeX*}  where  X*£  X  is  a 
finite  set  of  flows.  Then  the  |X*|  Okamoto  confidence  intervals 
that  result  from  (14a)  hold  simultaneously  with  probability  of  at 
least  a  when  1  —  (  1— a)  /  |X*|  replaces  a  in  (14a).  Also, 
confidence  intervals  based  on  normality  (14b)  hold  asymptotically 
when  this  same  substitution  for  a  is  made  in  (15).  These 
simultaneous  interval  results  follow  from  a  Bonferroni  inequality. 
See  Miller  (1981,  p .  8)  . 

3 .  Characterizing  Deterioration  Through  Time 

At  least  two  distinct  characterizations  exist  for 
multi-state  systems  with  random  capacities  for  which  maximum  flow 
is  a  performance  measure.  One  representation  views  system 
components  as  experiencing  successive  moderate-to-small  capacity 
reductions  of  random  magnitude  at  a  sequence  of  random  points  in 
time.  The  second  views  the  components  as  either  operating  at 
full  capacity  or  as  failed  at  random  points  in  time.  The  present 
paper  focuses  on  the  first  characterization  and  adopts  a 
formulation  from  the  theory  of  cumulative  processes,  as  described 
in  Smith  (1955),  to  facilitate  the  estimation  of  {l(x,t)}. 

Let  Tij  denote  the  time  at  which  the  ith  capacity 

reduction  occurs  to  component  j  so  that  t  ^  j  <  i|<j  k  =  i  +  1,  i+2 . 

Assume  that  (x.  -t.  .  Tn.=0;  i=1,2,...}  j=1,...,q  form  q 

independent  stationary  renewal  processes  with 


var 


<TiJ~Ti-1 


<  ® 


We  conceive  of  as  the  time  of  the  ith  degradation  for  arc  i 

in  arc  set  Ej  where 

Ej  n  Ek  =  0  j  *  k 
E  =  E.+. . ,+E  . 

1  q 

Define 


=  reduction  in  capacity  on  arc  j  at 
time  T^j 

and  for  each  j  take  A i j ,  to  be  independent  and 

identically  distributed  (i.i.d.)  random  variables  with 

EAij  =  Yj  <  - 

var  A ^ j  =  Xj  <  ® 

=  orr(»u>Alk)  -  »Jk  if  jEEk 

=  0  otherwise 


and 


cor r ( A . 

i  J 


T  i  k  Ti- 1  , k  } 


“j  if  jeEk 
0  otherwi se  . 


Let 


J(x,Ek)  »  number  of  degradation  times  in  [0,x]  in  set  Ek 
so  that  at  time  t  the  accumulated  degradation  on  arc  j  in  set 


E 


J(t  ,E  ) 

Y . ( t )  =  l  A.  if  J(t,E  )  >  0 

J  1=1  XJ  K 

=  0  otherwise 

and  the  capacity  of  arc  j  is 

B  (t)  =  B  (0)  -  Y j ( t ) . 

As  t  -»  « ,  one  has  (Smith  1  955  )  for  arcs  j.meE^ 

Y  ( t )  =  m  .  x  +  o(x), 

J  J 

where  o(t)  denotes  a  function  h  (  t  )  such  that  lim  h(t)/x  =  0, 

f  -►  CO 

"j  ■  Vv 

var  Y  ( t  )  =  o .  . t  +  o(t), 

sJ  J  \J 

°jj  =  Cxj  "  2wj(Yj/0k)(Bkxj)4  +  Bk(Vak)2]/ak 

and 

cov  [Y  (t),Y  (t)]  =  o,  t  *  o ( x )  if  j  e  E 
J  m  j  m  m 

=  0  otherwise 


where 


=  Cp  -i_t  A  jA  )2 
J  m  j  m 


“jtj  tBkAj ' 


Moreover  , 


a  s 


the  joint  distribution  of 


|[B.(t)  -  B.{0)  +  p.x]/x2  jeE}  converges  to  the  multinormal 
J  J  J 

distribution  with  mean  vector  0  and  covariance  matrix 
E  =|)  o.j  |  | .  Then,  for  large  t  { Z ( t , C )  Cerj  are  approximately 
normal  with 

EZ(t,C)  =  l  [B  (0)  -  p  t] 
jeC  J  J 

and  C , D e  r  • 

cov  [ Z ( t  ,  C  )  ,  Z ( x , D  )  ]  =  x  l  o., 

i , j  eC nD  J 

When  this  formulation  in  terms  of  cumulative  processes 
applies,  it  offers  many  conveniences  for  Monte  Carlo  sampling. 
Observe  that,  because  of  the  asymptotic  normality,  the  sampling 
distribution  of  B i ( x ),...,  Bn  (  x )  requires  knowledge  only  of  the 
initial  capacities  {  B  ^  (  0  )  1<i<n},  the  means  1  <  i  <  n  }  and 

the  covariances  l°ij  1£i,j£n}.  In  particular,  note  that  the 
covariances  i^ij  i * j  }  completely  characterize  statistical 
dependence  in  capacity  changes  on  different  arcs.  It  is  entirely 
plausible  that  one  can  derive  suitable  empirical  estimates  of 
these  quantities  from  time  series  or  cross-sectional  data  for 
multi-state  systems  encountered  in  practice.  Moreover,  sampling 
I Bj ( x )  j  e  C  }  subject  to  H(x,C,x)  is  easily  realized  using  the 
conditional  properties  of  the  normal  distribution  and  Algorithm 
RSNB  in  Cheng  (  1  9  8 1  ,  pp  .  3  1  M  -  3  1  5  )  ,  which  takes  0(|c|)  time. 

Two  restrictions  need  to  be  kept  in  mind  when  assessing  the 
suitability  of  this  probabilistic  model.  They  are 
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max  [B.(t)-B.(0)]<0  (18) 

jeE  J  J 

and 

min  B . ( x)  >  0,  (19) 

jeE  J 

which  reality  demands  with  probability  one  but  to  which  normal 
theory  assigns  some  positive  probability  less  than  unity.  In  the 
case  of  statistically  independent  capacities  this  probability  is 
1  -  0 ( t )  where 

n 

0( t )  =  1  -  n 
j  =  1 

and  <P  (  •  )  denotes  the  standard  normal  distribution  function. 
Clearly  1  -  0(x)  should  be  small  before  one  accepts  the 

suitability  of  the  normal  approximation.  As  an  additional 
precaution  when  using  the  normal  model  one  can  use 
max  1 0,  min[ Bj ( 0  )  ,  B  j  (  t  )  ]  }  in  place  of  Bj(t)  for  the  capacity 
of  arc  j  at  time  x  in  Procedure  MFBOUNDS. 

A .  Credibility  Analysis 

There  remains  the  question  of  how  large  K  should  be  in 
order  to  ensure  the  level  of  statistical  accuracy  required  for  a 
particular  study.  As  examples,  one  may  want  to  estimate  L(x,t) 
subject  to  the  absolute  error  criterion 

1 1  (  x  ,  e  )  =  {|  Lk  (x  ,  t  )-L  (  x  ,  x  )|  <e  }  e  >  0  (21) 


or  the  relative  error  criterion 


-23- 


r] 


I  2  (  x  ,  w  )  =  (|  Lk(x  ,  T  )-L(x  ,  T  )|  <0>L(x  ,  T  )  J  0<  u»<  1  .  (22) 

Some  studies  may  demand  that  the  least  criterion 

I-(x,e,u)  =  I  1  ( x  ,  e  )  u  I  2 ( x , u  )  (23) 

be  satisfied  and  others  may  demand  that  both  criteria 

Ilt(x,G,t*>)  =  I  1  (  x  ,  £  )  n  I  2  (  x  ,  w  )  (  2  <4  ) 

be  met.  Here  I  3  (  x  ,  e  ,  w  )  is  the  least  stringent  and  I  n  (  x  ,  e  ,  co  )  is 
the  most  stringent. 

Since  S(x)  in  Procedure  MFBOUNDS  has  the  binomial 
distribution  with  parameters  K  and 
p(x)  =  L (x  ,  r  )/N ( x  ,  t  ) , 

satisfying  I ]  ( x , c  ,u  )  =  I^tx.e)  is  equivalent  to  requiring  that 

* 

u ( x > e [ p  ( x  )  ,  u  (x)  ],  where 


U 1  * ( x  )  =  [S(x)/K-e/N(x,T  )  ] 


U1 ( X )  =  1  - [ 1 -S( x  )/K-e/N( x  ,  T  )  ]  , 

and  satisfying  I  2 ( x , e , u )  =  I  2  (  x  ,  w  )  is  equivalent  to  requiri 

* 

that  u ( x  )  e [ p  ( x  )  ,  p  ( x )  ]  where 


u2„(x)  =  S ( x ) / K ( 1 +u) 


vv 


T'  v’v 
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U*(x)  =  1— [1— S(x)/K(1— tu)  ]  +  . 


it 

Moreover,  satisfying  I  .  (  x  ,  e  ,  u  )  implies  p  (  x  )  e  [  u  ^  #  (  x  )  ,  p  ^  (  x )  ]  ,  wliere 


u3„(x)  =  min[Wl #(x)  ,u2*(x)] 


and 


(27) 


*  *  * 
u  (x)  =  max[ p 1 ( x  ) , u2( x )  ] , 


it 

and  satisfying  iJx,  e  ,01)  implies  p  (  x  )  e  [  p  ^  #  (  x  )  ,  p  ^  (  x  )  ] 


where 


l*  |j  *  (  x  ) 


=  max[p1 ^(x)  ,p2,(x)  ] 


and 


(28) 


* 

Uj,(  x ) 


*  * 

min[u1 (x) ,  p2( x ) ] . 


That  is,  for  S=S(x) 

1 -pr[I  .  (x  ,e  ,u>)  ]  =  pr {p( x)£[p  .  *( x)  ,p  .  ( x)  ]  } 

J  J  J 


(29) 


=  1-  F_(K,p.,(x))+F_(K,p.  (x))  1SJSH, 


S'  ’ 


S'  ’  H j 


F  being  the  binomial  distribution  function. 

For  sets  of  flows  X  =  {x1<...<xr}t  absolute  errors 

(e^i . cr}  and  relative  errors  { w wr } , 


1  -  pr  [  fT  I  (  x  ,  e  .  ,U).  )  ]  £  A, 

i  =  ]  J  1  1  1  J 


I 


>1 


i 


« 

-T-i 


;-.i 


where 


1.  -  l  [1~F_(K,v.<(x  ) ) +F„(K , y  (x. ))] 

J  i  =  1  5  j  1  b  J  1 


1 Sk£4  , 


the  upper  bound  following  again  from  a  Bonferroni  inequality 
(Tong  1980,  P.1H3). 

To  guide  the  choice  of  sample  size,  one  can  increase  K  in 
steps,  compute  the  bound  corresponding  to  the  chosen  criterion  at 
the  end  of  each  step  and  use  the  sequence  of  such  values  as  a 
measure  of  the  credibility  of  having  satisfied  the  criterion.  To 
ensure  some  progress  in  the  value  of  the  bound,  increments  in  K 
should  be  relative.  Doubling  K  at  each  step  is  a  reasonable 
choice.  IMSL  (1982)  and  SAS  (1982)  provide  routines  for 
evaluating  the  binomial  distribution  with  high  accuracy.  Note 
that  in  contrast  to  the  discussion  on  confidence  intervals  in 
Section  3,  the  development  here  relies  on  the  exact  binomial 
distribution. 


5 .  Estimating  Variance  Reduction 

As  Section  3  notes,  (9)  indicates  the  benefits  of  variance 
reduction.  One  can  estimate  this  quantity  by 


R  (  x  ,  t  )  =  -  x  Q  (  x  ,  t  ) 

T 

where 


Q(x,t)  =  [1-L„(x,t)]/[N(x,T)-L„(x,-t)], 


T  and  T  being  the  sample  mean  times  per  replication  for  Procedures 
MFCRUDE  and  MFBOUNDS  respectively.  Although  T/T  generally  is  a 


relatively  stable  quantity,  Q(x,i)  may  not  be,  especially  when 
L ( x  ,  t  )  is  close  to  N(x,t).  Therefore,  relying  solely  on  (26) 
as  an  indicator  of  variance  reduction  can  be  misleading.  To 
resolve  this  problem,  one  can  compute  a  confidence  interval  for 
R ( x , t )  in  (9) . 

Let  N  =  N ( x  ,  t )  ,  L  =  Lk(x,t)  and  R  =  R(x,t),  so  that 


Y  =  1-L  -  R(N-L) 


has  mean  zero  and 


var  Y  =  ( NR  -  1 ) ( 1~N)/K 
Then  Chebyshev's  inequality  gives 


pr 


|  l-L-R(N-L)  |  1 

C  (  NR  -  1  )  (1-N)/K]1/z  ( 1  -  a ) 


>  a 


yielding  a  100  «  a  confidence  interval 

(  1  -  L  )  (  N  -  L  )  +  B2N(1-N)/2K  ±  S  (  1  ~N  )  [  (N-L)  L/K  + 

(N-L)2 


2  2  2  V 

B  NV4K  ]2 


(27) 


i/ 

with  B  -  1/(1-a)  2.  If  warranted,  a  tighter  normal  interval  with  B 
computed  as  in  (15)  can  be  used. 


6 .  Example 

Figure  1  shows  a  directed  network  for  which  we  wish  to 
estimate  { L  ( x  ,  t  )  }  for  t  =  500  and  X  =  {800  +  1  0  0  i  i  =  1 . 12). 
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i 


Table  1  gives  the  corresponding  arc  parameters.  Since 


Insert  Fig.  1  and  Tables  1,  2  and  3  about  here. 


■»  -  »  1 


p 


s 

■ 


0(500)  =  .0172  in  (20)  the  normal  model  is  reasonable.  The 

selected  set  of  edge-disjoint  cutsets  is  C  •)  =  {  1  ,  2  .  3  }  » 

C2  =  { 1 ^ , 2H , 25  }  ,  and  C3  =  (5,6,9,10,12}.  Table  2  shows  (c(x.,t), 
N  (  x  ,  t )  1£i£r}.  Table  3  shows  the  results  of  using  Procedure 
MFB0UNDS  for  K  =  65536  independent  trials.  In  addition  to  the 
favorable  variance  reduction  for  all  x,  note  the  wide  confidence 
intervals  for  R  (  x  ,  t  )  for  900  £  x  £  1  400 ,  reflecting  the 
closeness  of  the  upper  bound  N ( x  ,  t  )  to  L  (  x  ,  t  )  .  Also,  note  the 
substantial  variance  reduction  for  x  £  1800  and  the  fact  that 
variance  reduction  is  not  monotone  in  x. 

Figure  2  provides  an  understanding  of  this  nonmonotone 
behavior.  These  cutsets  act  as  the  principal  restrictions  on 


Insert  Fig.  2  about  here. 


capacity  when  looking  at  flows  in  the  interval  [0,  2000].  For 
small  x,  Ci  is  the  dominant  restriction.  As  x  increases  C-|  and 
C  2  become  equally  dominant  and,  as  x  continues  to  increase,  C2 
becomes  the  restrictive  cutset.  The  point  of  lowest  variance 
reduction  apparently  occurs  where  Cj  and  C2  are  approximately 
equally  dominant. 


% 


Table  4  gives  individual  Okamoto  and  normal  .95  confidence 
intervals  for  {L(xitx)  1£i£r)  along  with  the  sample  skewness 
based  on  (16).  Although  both  sets  of  intervals  give  at  least  two 
digit  accuracy,  the  sample  skewness  encourages  one  to  rely  more 
on  the  Okamoto  bounds  for  900  £  x  £  1300. 


Insert  Tables  4  and  5  about  here. 


To  demonstrate  the  credibility  analysis  of  Section  4,  Table 
5  lists  the  bounds  for  the  four  criteria  described  there  for 
Ej  =  c  =  .01  and  wj  =  w  =  . 1  1£i£r  for  Procedure  MFBOUNDS  and  for 
crude  Monte  Carlo  sampling  in  Procedure  MFCRUDE.  Most  notable  is 
the  relative  error  criterion  which  shows  that  a  bound  less  than 
.1  is  achieved  with  high  certainty  at  K  =  4096  in  contrast  to 
Procedure  MFCRUDE  which  does  not  achieve  this  bound  with  moderate 
assurance  even  for  K  =  1048576.  That  is,  MFCRUDE  requires  a 
sample  size  at  least  102 4  times  larger  than  that  for  MFCRUDE  to 
achieve  the  same  relative  error  criterion. 

7 .  Essential  Steps  for  Implementation 

Here  we  list  in  appropriate  order  the  essential  steps  needed 
to  implement  the  proposals  in  this  paper; 

1.  Select  the  time  t  of  interest. 

2.  Select  the  flow  set  X  =  { x q  =  0< x i <  .  .  .  < xr  }  . 


-,1.  ->  -  \r*.™7%?K'"'e.''y7'\!ry,',.'+r!!V'.  >J  ->,  “''’VJ1’.1  1  *  ».«  "  »■  T7  ».» ■.»  V  r— T 
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3.  Compute  0(t)  in  (20)  to  check  the  appropriatenesss  of 

the  normal  model.  If  inappropriate,  stop. 

H .  Select  an  error  criterion  and  specify  the  corresponding 
{  e  . . er)  and  {  uj  -| . wr  }  (Section  *0  . 

5.  Determine  a  set  of  edge-disjoint  cutsets 

r  =  1  c  i . c  j  1 . 

6.  For  each  xeX,  determine  C(x,t)  and  N(x,t)  using 

Procedure  A. 

7.  Use  Procedure  MFBOUNDS  for  the  sampling  experiment. 

8.  If  the  experiment  is  to  be  performed  in  blocks  of 

K  -j  ,  K 1  +  K2,  Ki  +  K2  +  K 3...  r  epl  i  ca  t  i  ons  ,  then  after 

each  block  compute  the  credibility  probability  bound 
for  the  selected  error  criterion.  Suggested  increments 
are  Kj_  =  K  ^  2i_ 1  for  i  =  1  , 2  ,  .  .  . 

9.  After  completion  of  the  experiment  compute  the  sample 
skewness  in  (16)  and  a  confidence  interval  for  each 
L(x,t)  xeX  x*xQ. 
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Table  1 

Capacity  Parameters  for 
Network  in  Figure  1 


868.33 
2532.62 
1 736.66 
3473. 31 
361 8.03 
723.61 
3256.23 
21  70. 82 
1 736.66 
723.61 

868.33 
3256.23 
1 302.49 
1 809.02 

863.33 
2604.98 
361 8.03 
2894 . 43 
1 085. 41 
1 447.21 
3256.23 
2532 . 62 
1 809 . 02 
1 736.66 
21  70 . 82 


20 
50 
40 
80 
00 
00 
50 
00 
40 
00 
1  .  20 
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1  .  80 

2.50 
1  .  20 


I'I'J 


1  .5 

2.0 

4.5 

3.5 
2.50 
2. 40 
3.00 


14.40 

122.50 
57.60 

230.00 

250.00 

10.00 

202.50 
90. 00 
57.60 
10.00 

14.40 
202.50 

32.40 
6.25 

14.40 


129.60 
250.00 
160.00 
22.50 
4.00 
202 . 50 
122.50 

62 . 5 

57.6 
9.0 


Table  2 


Relevant  Disjoint  Cutsets* 
(t  =  500) 


i 

xi 

C  (  X  i  ,  T  ) 

N( Xi  ,  i  ) 

1 

900 

2.3} 

.  98627 

2 

1  000 

.  97024 

3 

1  1  00 

. 941 05 

4 

1  200 

. 89306 

5 

i  300 

.  821 81 

6 

1  400 

.  72628 

7 

1  500 

.61 062 

8 

1  600 

(14,20,25} 

.47865 

9 

1  700 

1 1 

2.3} 

.  281 63 

1  0 

1  800 

1  1  4 , 2( 

3,25} 

.03775 

1  1 

1  900 

. 0041 4 

1  2 

2000 

1 

.  0001  8 

.  T  )  } 

is  computed  as 

in  (11).  1 N ( x i , t ) , 

[  is  computed 

(12). 


flow 

upper 

bound 

estimate 

estimated 

variance 

.95  confidence 
interval 
for  R(x,r) 

X 

N(x.t) 

Lk(x,t) 

V(Lk(x,t)) 

R(x,t  ) 

T 

—  x  r(x,t)  lower  upper 

T 

900 

.9863 

.9861 

.31  70 

X 

10"8 

66.16 

22.21 

21  .98 

203.43 

1000 

.9702 

.9696 

.9856 

X 

10"8 

45.68 

15.33 

24.22 

86.99 

1 1  00 

.941  0 

.9399 

.1668 

X 

10-7 

51  .69 

17.35 

32.01 

83.88 

1200 

.8931 

.8906 

•  3370 

X 

10"  7 

44.12 

14.81 

32.02 

60.96 

1  300 

.8218 

.8172 

.5754 

X 

10"  7 

39.61 

13.30 

31.63 

49.70 

1  400 

.7263 

.7177 

.9358 

X 

10"  7 

33.04 

11.09 

28.30 

38.60 

1  500 

.6106 

.5780 

.2879 

X 

10-6 

12.93 

4.34 

12.08 

13-84 

1600 

.4787 

.3670 

.6252 

X 

10“6 

5.67 

1.90 

5.52 

5.82 

1  700 

.2816 

.1392 

.3025 

X 

10"6 

6.04 

2.03 

5.96 

6.13 

1800 

.3775  * 

10'1 

.2411  x 

10"1 

.5016 

X 

10“8 

71  .60 

24.03 

69.98 

73.26 

1  900 

.41 42  * 

1  O'2 

.1626  * 

10“2 

.6243 

X 

10" 10 

397.02 

133.27 

391.56 

402.69 

2000 

.1788  * 

10“  3 

.3900  x 

io-1* 

.8321 

X 

1  0“ 1  3 

7151 .79 

2400.65 

7087.20 

7219.18 

T  and  T  were  obtained  from  runs  with  MFCRUDE  and  MFBOUNDS  respectively. 
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